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We report on a first relativistic calculation of the quenched bottomonium spectrum from 
anisotropic lattices. Using a very fine discretisation in the temporal direction we were able to 
go beyond the non-relativistic approximation and perform a continuum extrapolation of our results 
from five different lattice spacings (0.04 — 0.17 fm) and two anisotropies (4 and 5). We investigate 
several systematic errors within the quenched approximation and compare our results with those 
from non-relativistic simulations. 
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The recent experimental activities at B-factories have triggered many theoretical attempts to understand heavy 
quark phenomenology from first principles. The non-perturbative study of heavy quark systems is complicated by 
the large separation of momentum scales which are difficult to accommodate on conventional isotropic lattices. These 
problems are particularly severe in heavy quarkonia, where the small quark velocity separates the high momentum 
scales around the heavy quark mass from their small kinetic energy; m q 3> m q v 2 . While this motivates the adiabatic 
approximation and potential models, a full dynamical treatment poses significant computational problems. The 
\ limitations of present computer resources force us to tackle this problem in some modified framework. Thankfully 
0^ ■ there is also a wealth of spectroscopic data against which different lattice methodologies and improvement schemes 
can be tested. Once experimental results can be understood with some accuracy the lattice will provide a powerful 
tool for other non-perturbative predictions such as decay rates and form factors involving heavy quarks. 

Several approximations to relativistic QCD have been proposed to describe accurately the low energetic phenomcnol- 
I , ogy of heavy quarkonia However, the (non-)perturbative control of systematic errors in those approximations 

C^) • is very difficult and in practise one has to rely on additional approximations. The high precision results for the spin 
structure in non-relativistic bottomonium calculations || have been hampered by large systematic errors which are 
C^j i difficult to control within this effective theory. Higher order radiative and relativistic corrections are sizeable for 
bottomonium [Q|5| and even more so for charmonium Q . Still more cumbersome are large observed scaling violations 
Oh, HQ that cannot be controlled by taking the continuum limit. 

We take this as our motivation to study heavy quarkonia on anisotropic lattices in a fully relativistic framework. 
This approach has previously been used to calculate the charmonium spectrum with unprecedented accuracy ^ [LlJ . 
Here we extent the method to even heavier quarks and focus on our new results for bottomonium. Preliminary results 
k> | were already reported in [|l0| . 

There are additional advantages in using anisotropic lattices which have been employed previously. It has been 
demonstrated by several authors that a fine resolution in the temporal direction is useful in controlling and extracting 
higher excitations and momenta in the QCD spectrum [fl2| fl5[ . It has been noticed long ago that anisotropic lattices 
are also the natural framework to study QCD at finite temperature This has been used in |[7|,[l8) to retain 

more Matsubara frequencies at high temperature. Furthermore, it has been suggested that anisotropic lattices can 
circumvent problems arising due to unphysical branches of the dispersion relation on highly improved lattices [ 1^] . 

In Section [il] we introduce the anisotropic gauge and quark action for our study. In Section |l| we give the details 
of our simulation and study systematic errors in Section IV . Our results are discussed in Section [v| and Section VI 
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concludes this paper. 



II. ANISOTROPIC LATTICE ACTION 



The basic idea of our approach is to control large lattice spacing artifacts from the heavy quark mass by adjusting 
the temporal lattice spacing, at, such that m q a t < 1. For the bottomonium system this implies a^ 1 > 5 GeV. Dis- 
cretisation errors from spatial momenta are controlled by \p\a s < 1, which can be satisfied more easily at conventional 
spatial lattice spacings aj 1 > 1.5 GeV. We retain the isotropy in all the spatial directions. The continuum limit can 
then be taken at fixed anisotropy, £ = a s /at- 
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More specifically, we employ an anisotropic gluon action, which is accurate up to 0(a 2 , a 2 ) discretisation errors: 



s = -P ( E ^P^) +^oPit(x) I . (1) 

l a;,i>j x,i 

This is the standard Wilson action written in terms of simple plaquettes, P fiv {x). Here £o) are two bare parameters, 
which determine the spatial lattice spacing, a s , and the renormalised anisotropy, £, of the quenched lattice. For the 
heavy quark propagation in the gluon background we used the "anisotropic clover" formulation as first described in 
1^,^). The discretised form of the continuum Dirac operator, Q = m q + Ip, reads 

a s 

Q = mo + vt WqJq + v s Wiji - — [ct crokFok + c s (T k iF k i\ , 
W„ = - (o^/2) 7m A m . (2) 

Here denotes the symmetric lattice derivative and a 2 A M g(a;) = U fl (x) q(x + fi) — 2 q(x) + U- I _ l (x) q(x ~ fi). For 
the electromagnetic field tensor F^ v we choose the traceless cloverleaf definition which sums the 4 plaquettes centered 
at point x in the {p, v) plane: 

a^a-v F^(x) = - [Pp V (x) + P v p,{x) + Pp,p(x) + P PlM (x) - h.c] . (3) 

As has been demonstrated by P. Chen ||, Eq. (||) is equivalent to the continuum action (up to 0(a 2 ) errors) through 
0(3) symmetric field redefinition: 

q(x) -> ( 1 + ^-a t m q + ^a t f t + ^a a f^j q(x) . (4) 

In a recent paper po| , S. Aoki et al. observed that this 5-parameter fermion action does not generate Euclidean 
0(4) covariant quark propagators. While it is true that the fermion Green's functions derived from our action will 
contain non-covariant terms at 0(a), these can be removed by undoing the field transformation (0) at the end of the 
calculation. Rather than working with a more complicated action where additional parameters have to be tuned for 
covariance, it would be easier to implement this change of spinor basis after the expensive inversion of the fermion 
matrix. However, in this work we do not even need such a final transformation since we do not study spin-1/2 hadrons. 
[] We have chosen Wilson's combination, W M , of first and second derivative terms so as to ensure the full projection 
property and to remove all doublers. The five parameters in Eq. (^) are all related to the quark mass, m q , and the 
gauge coupling as they appear in the continuum action. By tuning them appropriately we can remove all 0(a) errors 
and re-establish Euclidean 0(4) symmetry. The classical estimates have been given in M: 



m = m q {l + - a s m q ) (5) 



1 + g a s m q 
1+2 a t m q 

c s = v s (7) 
ct = Uu s + uA) . (8) 



Simple field rescaling enables us to set one of the above coefficients at will. For convenience we fix v s = 1 and adjust 
v t non-perturbatively for the mesons to obey a relativistic dispersion relation ( c(0) = 1 ): 

£; 2 (p)=i? 2 (O)+ c 2 (p)p 2 + 0(p 4 )... . (9) 

We also chose mo non-perturbatively, such that the rest energy of the hadron matches its experimental value 
(M( 3 S7~) = 9.46 GeV for bottomonium). For the clover coefficients (c s ,ct) we take their classical estimates from 
Eqs. (ffl) and d§) and augment the action by tadpole improvement: 



1 We would like to thank N. Christ for clarifying this point 
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Ui -» U t /u s (10) 
U ^U /u t . (11) 



Following Ref. Q the resulting (dimensionless) fermion matrix can then be written as 



v s - 1 

a t Q = u t m + v t W 7o + t~ W iH ~ o 

SO ^ 



— ^ cr feF fe + -^r- a k [F k i I , (12) 



where we arranged the terms in such a way that the temporal Wilson term is multiplied by only v t . To compute the 
coefficients of the spatial terms we used £/£o = u t /u s , which is accurate within 3% [PJ?lf ■ Remember that we can 
still choose v s at will and a 3% uncertainty for c s is certainly acceptable - at our lattice spacings a different tadpole 
prescription has a much more pronounced effect on this coefficient. The tadpole coefficients have been determined 
from the average link in Landau gauge: = 1/3 (tr Uy,(x)) Landau- For brevity we will refer to this as the Landau 
scheme. Any other choice for (c s ,Ct) will have the same continuum limit, but with this prescription we expect only 
small 0(aa) discretisation errors. 

III. SIMULATION DETAILS 

For the generation of quenched gauge field configurations we employ a standard heat-bath algorithm as it is also 
used for isotropic lattices. The renormalised anisotropy, £, is related to the bare parameter £o through: 

£ = »?(&/3Ho • (13) 
A convenient parametrisation for r) can be motivated by a one-loop analysis of the renormalised anisotropy: 

£ (3 1 - a g 



^/?) = 1 + (1-7)^T-^ , (14) 



where ?7i(£) ~ 1.02 and almost constant for our range of £. We use Pade parameters (ao = —0.77810 and a% — 
—0.55055) which have been determined by an excellent fit over a wide range of lattice data with different /3 and £ 
plf . Based on Eq.|l4| we choose the appropriate £o to realize £ close to 4 and 5 for each value of f3 = 6/g 2 — 5.7 — 6.5 . 
Depending on the gauge coupling, we measure hadron propagators every 100-400 sweeps in the update process, which 
is sufficiently long for the lattices to decorrelate. 

The construction of meson operators with given J PC assignment is standard. As in Rcfs. ]l0| , p2[ , we use both local 
and extended operators by employing the 16 T-matrices and spatial lattice derivatives for the quark-bilinears: 

M{x)=q{x)T i A j [U\& k [U}q{x) . (15) 

To improve the overlap with the ground state we also implemented a combination of various iterative smearing 
prescriptions for the quark fields and gauge links: 



(x) -> q {n) (x) = (1 - -A 2 [U]) n q(x) Jacobi smearing (16) 



U{x) -> U {m+1 \x) =V SU (3) \U {m) {x) + a ^ S(x)\ APE smearing (17) 

staples 

For conventional states, we find it often sufficient to use (gauge-fixed) box sources of different size. We use local 
sinks throughout. This setup allows us to extract reliably both the ground state energies and their excitations from 
correlated multi-state fits to several smeared correlators, C s (t), with the same J PC : 

C s (t, p) = (M(t, p) M s (0, p)> = J a|(p) (e-^ (p) 4 + e-^W"'-*) . (18) 

i=l 

Since we are working in a relativistic setting, the second term takes into account the backward propagating piece from 
the temporal boundary. An example of this is given in Figure [l], where we show the effective masses for 3 S± = 
from three different box sources along with the fit results from a 2-state and 3-state ansatz. In this example the 
smallest box shows clear contributions from excited states, while the other two sources project much better onto the 
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ground state. We did not further optimize the box size for perfect overlap, but since for the largest source the plateau 
is approached from below we can estimate that an "optimal" box would have a physical extent of 0.2-0.3 fm. This is 
in good agreement with phenomenological expectations about the size of the T ground state. What is more important 
here for us is to have a variety of different overlaps which are well suited to constrain the multi-state fits. 

To obtain the dispersion relation, E(p), we project the meson correlator also onto four different non-zero momenta 
by inserting the appropriate phase factors, exp (— ipx), at the sink. This is illustrated in Figure where we use 

p=0£)n with n= (0,0,0), (1,0,0), (1,1,0), (2, 0,0), (2,2,0) . (19) 

For an estimate of the spin-splittings we take the difference of the ground state masses from multi-state fits. Because 
these states are highly correlated we obtain very accurate estimates for the spin-splittings compared to the absolute 
masses. An example of this is shown in Figure^. We use jackknife and bootstrap ensembles for our error estimates and 
find good agreement in all cases. In order to call a fit acceptable we require consistency within errors for different fit 
ranges and Q- values to be bigger than 0.1. Our main results are collected in Table | along with the input parameters 
for the different data sets. 

For (/?,£) = (6.1,4) and (6.3,4) we have also analyzed extended operators that give access to 2 ++ = 3 p2, -D-states 
with two units of orbital angular momentum and exotic hybrid states, Table [n| These are more states than previously 
obtained from non-relativistic calculations. Even though our scaling window is not wide enough to make a serious 
attempt of a continuum extrapolation, we expect lattice artifacts to be less pronounced for those large states. At 
(6.1,4) we find a good agreement of our results and those reported in ]l3| , |l4| , p2"| . This is not unexpected since the 
non-relativistic approximation is known to work better for hybrid states, which are thought to live in a very shallow 
potential [H. In particular we find the lowest lying hybrid to be the 1 '". We are however concerned that these 



results are more sensitive to finite volume effects than those for S and P states, see Sec. [V. Especially at (6.3,4) we 
have reason to believe that our lattice is not sufficiently big to accommodate such large wavefunctions, resulting in an 
inversion of the characteristic level ordering. A simulation on bigger volumes is necessary to simulate hybrid states 
more accurately. 

In the following we will focus on the low-lying bottomonium states and their spin-structure where statistical errors 
are small enough to allow a critical comparison with results from the non-relativistic approximation. 



IV. SYSTEMATIC ERRORS 



Apart from our main data sets we performed extensive checks of systematic errors within the quenched approx- 
imation. As mentioned above, we enforce a relativistic dispersion relation by tuning the bare parameter v t non- 
perturbatively at each quark mass to obtain c=1.0 very accurately within 1 — 2%. From a classical analysis of the 
dispersion relation we expect c oc 1/vt- Indeed, for an arbitr ary decrease of v t by 10% we observe an increase of the 
velocity of light by a similar amount (7.8(2.1)%), see Table III and Figure ^. Within the statistical error for the P 
structure we could not resolve any significant effect, while the hyperfine splitting showed an increase by 8(10)%. From 
those observations, we expect only small changes (< 2%, < 1 MeV) due to c ^ 1. 

We also tested for finite size effects by explicitly comparing the results from two different spatial volumes of (1.3 fm) 3 
and (0.7 fm) 3 at (/?,£) = (6.1,4), see Figure [|. While changes in the ground state masses amount to only a fraction 
of a percent, the hyperfine splitting shows a slight increase of 6(4)% and the 1 P\ — 3 Si splitting is reduced by 10(9)% 
when going to the larger lattice. Therefore we expect the hyperfine splitting to be systematically underestimated 
when L < 1 fm. For the fine structure we could not resolve any shift within the statistical errors. For all but our 
finest lattice we have L w 1 fm which should be sufficiently big to accommodate the small bottomonium ground states 
in a quenched simulation. 

Another source of systematic errors is the tuning of the bare quark mass parameter, mo- Especially for the hyperfine 
splitting we observed a strong dependence on the ground state mass, M(mo), as shown in Figure O. On our coarser 
lattices we were able to tune the mass very accurately to within 10% of the physical bottomonium mass, but as the 
tuning becomes more and more expensive on our finer lattices, the deviation from bottomonium reaches 35(15)% at 
(/3,£) = (6.5,4). We use the potential model prediction, Ai^hfs oc 1/M, to rescale all measured AEhf s to the physical 
point, M = 9.46 GeV. Our results in Table | include these adjustments in both the central values and their errors. 
Since all spin splittings vanish in the static limit we also apply the same rescaling technique to the fine structure. 

Finally there will be systematic effects at finite lattice spacing, which are absent in the continuum limit. The tree- 
level estimates for c s and Ct will experience a shift due to radiative corrections, thereby changing the magnitude of the 
0(aa) artifacts. In Figure ^ and Table IV we compare the Landau scheme with another popular scheme, where the 
tadpole coefficients are estimated from the plaquette values: 3u M = (tr -P^) 1 / 4 . This amounts to a (10%, 2%) change 
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in (c s ,ct) at (5.9, 4). The corresponding change in the hyperfine splitting is -15(17)%. While this may indicate some 
additional reduction at coarse lattices (resulting in a slightly worse scaling behaviour) , it is also clear that within our 
errors we are not very sensitive to this choice of tadpole coefficients. Since any remaining difference will vanish in the 
continuum limit we did not investigate this dependence further. 



V. DISCUSSION 



To convert lattice data into dimensionful quantities we used the lattice spacing determined from the Pi — S% 
splitting which is fixed to match the experimental value (« 440 MeV). It is well-known that, without dynamical 
sea quarks in the gluon background, the definition of the lattice spacing is ambiguous and one cannot reproduce all 
experimental splittings simultaneously. Non-relativistic lattice calculations with two dynamical flavours |j[7j resulted 
in shifts of up to 5 MeV for the hyperfine splitting, but they are also not completely free from ambiguities in the 
lattice spacing and from systematic errors such as lattice spacing artifacts, radiative and relativistic corrections. For 
the purpose of this paper we accept the shortcomings of the quenched approximation and aim to control the other 
systematic errors instead, the combined effect of which could also be as large as 20% for bottomonium on the lattices 
in this study. 

The main features of the bottomonium spectrum from relativistic lattice QCD are summarised in Figure |8[ where 
we can see clearly the characteristic level ordering as additional angular momentum is inserted into the ground state. 
The overall agreement with experimental data is impressive, as is the spin-structure predicted from first principles. 
While it is probably too early to investigate the spin-splittings in excited states (2 3 5i —2 1 So has still large statistical 
errors), we have rather accurate data for the spin-structure of the ground states. 

In Figure ^ we plot the hyperfine splitting against the spatial lattice spacing at fixed anisotropy and compare it 
with previous simulations in this region. Already at finite lattice spacing we can see significant deviation from non- 
relativistic simulations. There could be several reasons for this. First, we should expect different discretisation errors 
from the different gluonic actions that have been employed in the past. While the unquenched results in Figure ^ 
were obtained from an RG- improved gluon action with two dynamical flavours |7) , many other quenched simulations 
where done using the standard plaquette action [||-| 23 j. Non-relativistic results from coarse and anisotropic lattices 
also exist |24j | and they seem to show scaling violations similar to those from unimproved gluon actions, even though a 
plaquette-plus-rectangle (Symanzik) action had been used. This could be an indication of comparatively large 0(aa) 
corrections which can compete with 0(a 2 ). 

Some of those results mentioned also use a different prescription for the tadpole coefficients than in our present 
work. We want to stress again that the tadpole prescription is merely an empirical way of keeping radiative corrections 
small. It has been suggested by several authors that the Landau scheme is particularly suitable to account for radiative 
corrections due to tadpole diagrams ]q,p5|. In the continuum limit (U^ 1) both definitions (Landau and plaquette) 
will be identical, but at finite lattice spacing they result in different estimates for the clover coefficients. We have 
demonstrated above that a lower value for uq will result in a slightly larger hyperfine splitting. A more pronounced 
effect has been observed in In the non-relativistic framework one can also quantify this shift which has been 
employed in ^3|. These authors rescale their results to the Landau scheme and found improved scaling. Here we 
find some indication for a minor improvement, but none of the schemes is completely successful in removing the 
apparently large O(aa) errors. A non-perturbative determination of the clover coefficients is in order to reduce those 
errors significantly. We used the Landau scheme mainly for simplicity. 

In contrast, the fine-structure is not expected to be as much affected by the tadpole prescription which can explain 
the better agreement of our estimates for 3 Pi — 3 Pq with the NRQCD values in Figure |l0|. 

We interpret the remaining differences between our results and the quenched NRQCD values as due to higher order 
relativistic effects. It has already been shown in j|,[3| that those corrections can be sizable (> 10% for bottomonium at 
0{mv e )). However, it is not trivial to extend the non-relativistic approximation to ever higher order. Our calculation 
does not suffer from this uncertainty since we are working in a fully relativistic setting. 

Above all we are now able to perform a continuum extrapolation of the hyperfine splitting. As shown in Figure ^ 
the scaling violations are large and we parameterise our lattice results by 1) a linear ansatz, 2) by a quadratic ansatz 
and 3) by a linear-plus-quadratic fit in a s . All these fits are acceptable with Q = 0.3 . . . 0.9. For the final analysis 
we also decided to omit the results from our coarsest lattice where we should expect potentially large lattice spacing 
corrections from higher orders due to citm q > 1 and a s m q v > 1. Including this data, however, does not change our 
continuum estimates significantly. If we assume that our action successfully removes all 0{a) errors, we can quote 
51.1(3.1) MeV for the quenched hyperfine splitting in bottomonium. This may be too optimistic and, allowing for 
O(aa) errors, we find 59(20) MeV from a linear-plus-quadratic fit. As can already be seen from the large error on 
the intercept, such a general fit is not very well constrained by our data and even the quadratic term is consistent 
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with zero. For comparison, we find 58.7(5.5) MeV from a simple linear fit. Notice though, that all these fitting 
methods result in comparatively high values given phenomenological estimates of 30-40 MeV. It remains to be seen 
how unquenching will affect this quantity, although we expect that it would further increase the hyperfine splitting. 
An accurate experimental determination of the bottomonium hyperfine splitting would be most valuable to judge the 
reliability of our approach. 

For the fine structure, 3 Pi — 3 Po, the situation is very similar as shown in Figure Wc find a rather large value 
of 50.9(6.4) MeV from the simple linear extrapolation in a s (Q — 0.22). This is 3er above the experimental value of 
32.1(1.5) MeV. The quadratic fit is slightly worse (Q = 0.14) and predicts 41.7(4.1) MeV for the continuum limit. As 
for the hyperfine splitting, a fit involving both terms is badly constrained and gives 68(22) MeV with Q = 0.12. 

In Figure [ll] we also plot the higher excitation, 25 — 15, against the spatial lattice spacing. Naturally these 
results have larger errors than spin-splittings and ground states as they rely on the control of additional parameters 
(amplitudes and energies) in the fitting process. In line with other simulations, we find a continuum result which is 
slightly larger than the experimental value. However, it will be important to measure this quantity more accurately 
on bigger volumes as our lattices may not be big enough for those larger states. Based on the observed finite 
volume effects, Figure [| we have reason to believe that the result from our finest lattice (/3, £) = (6.5,5) could be 
overestimated. Therefore we treat our data for the 25 — 15 splitting with caution. 



VI. CONCLUSION 



In conclusion, we have demonstrated that also systems containing quarks as heavy as bottom can be treated 
relativistically within the framework of anisotropic lattice QCD. Imposing a fine temporal discretisation appears to 
be a very natural way to respect the physical scales in the problem and it should open up a reliable alternative to 
non-relativistic simulations (where systematic errors are much harder to control). In our study of the bottomonium 
system we found noticeable deviations of the hyperfine splitting from non-relativistic simulations, while the fine- 
structure agrees well with previous results. Overall we observed strong scaling violations resulting in a relatively large 
value and large errors for the continuum limit. This however, is not totally unexpected as we only used a tree-level 
estimate for the clover coefficients and did not tune them non-perturbatively to eliminate 0{a s ) effects completely. We 
believe that further studies on anisotropic lattices will be able to determine the continuum spectrum of bottomonium 
with increasing accuracy. A non-perturbative determination of the clover coefficients is highly desirable, but it remains 
to be seen whether this is feasible as suggested in pq |. Alternatively, simulations at even finer lattice spacings would 
be useful to disentangle remaining 0{a s ) and 0(a z s ) errors. What remains to be done is the efficient implementation 
of anisotropic sea quarks which is necessary to control the systematic errors from the quenched calculation. Work in 
this area is in progress |27]] . 

Future applications are not restricted to heavy quark masses alone, but they may also be extended to heavy-light 
systems where both quarks can now be treated within a uniform approach. 
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2.292804 


2.059794 


2.074632 


1.946351 


1.946351 


1.88063 


1.7680 


1.7054 


c t 


1.170549 


1.127661 


1.1177481 


1.098368 


1.098368 


1.05055 


1.0154 


0.9002 


c(0) 


1.044(25) 


1.026(32) 


0.990(11) 


0.979(7) 


1.018(14) 


0.994(11) 


1.004(19) 


0.992(20) 


3 Si [GeV] 


9.79(66) 


9.53(34) 


9.57(44) 


10.40(30) 


9.49(74) 


10.11(97) 


12.60(67) 


13.04(14) 


1 3 Si - x 5o [MeV] 


36.0(3.5) 


37.2(2.6) 


37.5(2.9) 


44.8(2.0) 


33.8(3.1) 


37.0(5.7) 


49.9(4.9) 


51.0(7.9) 


2 3 Si - 1 So [MeV] 


155(273) 


56(57) 


37(12) 


40(11) 


57(18) 


2(58) 


48(41) 


57(29) 


x Pi - 3 P [MeV] 


16.2(4.7) 


20.6(4.5) 


25.2(3.2) 


39.5(5.7) 


28.7(6.6) 


30.7(5.5) 


54.0(7.7) 


46.6(8.7) 


1 Pi- 3 Pi [MeV] 


-1.4(2.3) 


-1.3(1.1) 


4.6(1.3) 


8.0(3.3) 


7.9(4.9) 


3.7(2.0) 


8.8(3.5) 


8.4(3.8) 


3 Pi - 3 P [MeV] 


16.0(2.8) 


23.7(3.1) 


24.5(2.2) 


28.3(4.0) 


28.0(5.7) 


27.4(4.2) 


44.7(5.1) 


39.9(6.6) 


R2S 


1.00(23) 


1.37(18) 


1.40(14) 


1.27(14) 


1.38(17) 


1.30(25) 


1.53(18) 


1.53(20) 


Ras 


1.99(50)* 


2.35(29)* 


3.13(32)* 


2.27(34)* 


3.6(1.3)* 


3.05(57)* 


3.25(51)* 


5.17(76)* 


Rip 


2.21(48)* 


2.32(30) 


2.26(26) 


2.40(13) 


2.26(30) 


1.98(61) 


2.56(49) 


2.14(34) 



TABLE I. Simulation Parameters and Results. We label our runs by (/3, £), where /3 is chosen to vary the spatial lattice 
spacing and £ determines the appropriate bare anisotropy £o, see Eq. Il4l The temporal lattice spacing, at, is determined from 
1 Pi — 3 Si. Ground state and higher excitations are from correlated 3-cosh fits to all channels. For the excitations we quote their 
normalized results Rx = {X — IS) /(IP — IS). An asterisk denotes results which have not been checked by fitting even higher 
excitations than those listed. For example, we did not perform 4-cosh fits to obtain better estimates for 3S. All spin-splittings 
have been rescaled to the physical bottomonium ( 3 Si = 9.46 GeV). 



(AO 


(6.1, 4) 


(6.3 ,4) 


configs 


330 


360 


separation 


400 


400 


3 P> - 'Pi [MeV] 


24.4(4.1) 


20.6(6.9) 


3 P 2 - 3 Pi [MeV] 


39(14) 


27(11) 


2 — _ 2"+ [MeV] 




31.0(6.9) 


R 2 -+ 




1.63(11) 


R»— 


1.81(26) 


1.67(11) 


Ri-+ 


3.41(43) 


4.25(34) 


R +- 


4.41(85) 


3.37(72) 


R 2 +- 




5.93(45) 



TABLE II. Results from extended meson operators. For the two lattices shown all run parameters are exactly as in Table Q 
and Rx is the normalized excitation above the ground state: Rx = (X — IS) /(IP — IS). However, here we employ also the 
extended operators of Eq. ( (lB] ) which give us access to 3 P2, states with L = 2 and the exotic hybrid candidates (0 +_ , 1 '", 2 H ). 
For P-states (L = 1) we use a single symmetric lattice derivative to create extended operators and for D-states we take the 
anticommutator {Ai, Aj }, (L = 2). The commutator [A», Aj] amounts to an insertion of a colour-magnetic field into the quark 
bilinear and gives rise to hybrid states. Different Dirac matrices result in different J FC assignments. 
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{Vs,Vt) 
(Cs,C t ) 


(1, 1.815) 

(2.074632, 1.1177481) 


(1., 1.650) 

(2.074632, 1.090685 ) 


configs 


780 


500 


c(0) 

a' 1 [GeV] 


0.990(11) 
9.27(43) 


1.067(18) 
8.20(34) 


3 Si [GeV] 

3 Si - 'So [MeV] 

'Pi - 3 P [MeV] 

'Pi - 3 Pi [MeV] 

3 Pi - 3 P [MeV] 

R 2S 

R2P 


9.57(44) 

37.5(2.9) 

25.2(3.2) 

4.6(1.3) 

24.5(2.2) 

1.40(14) 

2.26(26) 


9.05(38) 

34.4(2.6) 

24.3(5.7) 

4.6(1.4) 

24.6(2.5) 

1.324(81) 

2.54(23) 



TABLE III. Change in vt. In this table we comapre our results from (/?, £) = (5.9, 5) when the optimal v t = 1.815 is changed 
arbitrarily by 10%. This entails a change in c t . All other run parameters are as in column 4 of Tab. |. The corresponding 
changes in the spectrum as listed above and the lattice numbers are given in Fig. m 



tadpole scheme 


Uol 


HOP 


uol 


Uqp 




(5.9,4) 


(5.9,4) 


(6.1, 4) 


(6.1, 4) 


(N s ,N t ) 


(8, 64) 


(8, 64) 


(16, 96) 


(16, 96) 


configs 


600 


600 


660 


130 


a t _1 [GeV] 


6.57(43) 


6.40(53) 


10.57(31) 


10.6(1.1) 


a a [fm] 


0.1202(80) 


0.123(10) 


0.0747(22) 


0.0745(77) 


£0 


3.139035 


3.139035 


3.210801 


3.210801 


UOs 


0.785945 


0.810698 


0.800927 


0.822785 


MOt 


0.986984 


0.949731 


0.988702 


0.952954 


a s m q 


1.1200 


1.1200 


0.6700 


0.877 


{Vs,V t ) 


(1,1.50) 


(1,1.50) 


(1,1.573) 


(1,1.573) 


Cs 


2.059794 


1.876818 


1.946351 


1.795316 


Ct 


1.127661 


1.101421 


1.098368 


1.079827 


c(0) 


1.014(15) 


1.023(14) 


0.979(7) 




3 Si [GeV] 


9.25(61) 


9.04(75) 


10.40(30) 


12.3(1.3) 


3 Si - 'So [MeV] 


35.1(3.9) 


30.6(3.6) 


44.8(2.0) 


43.7(6.6) 


'Pi - 3 P [MeV] 


26.0(3.9) 


25.0(4.0) 


39.5(5.7) 


41(14) 


'Pi - 3 Pi [MeV] 


1.3(1.7) 


1.7(2.3) 


8.0(3.3) 


9.6(8.2) 


3 Pi - 3 P [MeV] 


23.7(2.7) 


22.3(2.8) 


28.3(4.0) 


30.8(9.4) 


R2S 


1.45(13) 


1.41(13) 


1.27(14) 


1.37(23)* 


R2P 


2.22(49) 


1.81(28) 


2.40(13) 


3.14(97)* 



TABLE IV. Here we compare our results from different tadpole descriptions at two different lattices spacings. We distinguish 
between the Landau scheme, uol, as defined in Section [il] and the plaquette scheme, uop, for which 3it M = (tr Pi,,)'^ 4 . While 
there is some indication that the hyperfine splitting increases when uo is decreased at finite lattice spacing, there is no clear 
improvement of scaling. See discussion in main text. 
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FIG. 1. Effective masses and fit results. We show the effective masses, open circles, for 3 Si at f ) = (6.5, 4) from 3 different 
box sources with spatial extent of 3, 6 and 9 in lattice units. The two larger sources project better onto the ground state as 
expected. The results from the 2-state and 3-state fits are also shown as filled squares and diamonds, respectively. We observe 
stable fit results for large enough t m in- In this example, we fix t ma x = 80. 
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FIG. 2. Dispersion relation. On the left-hand side we show t m i n -plots for the fitted energy difference between non-zero 
momentum states and the Si at rest. On the right hand side we plot the fitted splitting against the spatial momenta squared, 
p 2 = (2tv/L) 2 n 2 . Since we tuned v t very carefully, this correspndis to a velocity of light c(0) = 0.990(11). This example is 
from (/?,£) = (5.9,5) with v t = 1.815. 
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FIG. 3. Here we show fit results for the 3 Pi - 3 P splitting at (/?,£) = (6.1,4), 16 3 x 96. We fix t max = 48 and vary i min 
in search for a plateau. Strong correlations between the the two states allow a very accurate determination of the energy 
difference. 
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FIG. 4. !/ t -tuning. We illustrate the non-perturbative tuning of the velocity light (/?, £) = (5.9,5). Our best estimate is 
v t = 1.815 and we show the effect on spectral qunatities as we change this parameter by 10%. 
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FIG. 5. Finite volume effects. Here we compare our spectrum results at (/3, £) = (6.1, 4), TVs x 96 for N s — 8 and 16. 
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FIG. 6. Mass dependence. We plot the dimensionless hyperfine splitting against the inverse of the measured meson mass. 
The vertical line denotes the Bottomonium mass at (/3,£) = (6.1, 4). The corresponding bare quark masses are mo = 0.67, 0.40 
and 0.10. 
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FIG. 7. Different tadpole prescriptions. Here we compare the effect of changing the value of (c 3 ,Ct) by (10%, 2%). Within 
our statistical errors there is only a minor shift to indicate the corresponding change in 0(aa) scaling violations. 
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FIG. 8. Relativistic bottomonium spectrum. We plot our quenched continuum estimates as full circles. Selected results for 
2 and exotic candidates (hybrids) from finite lattice spacing, (6.1,4), are shown as open squares. The 1 Pi — 3 Si splitting is 
used to set the scale. Where available, experimental values are shown as horizontal lines. 
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FIG. 9. Continuum extrapolation of the hyperfine splitting in Bottomonium. Circles denote our results from different 
anisotropic lattices. Open squares and triangles are quenched NRQCD results from isotropic and anisotropic lattices. Filled 
diamonds are unquenced NRQCD results with two dynamical flavours at m w /m p ~ 0.5. We also show two different continuum 
extrapolations of our £ = 4 results, where we assume a simple linear or quadratic form in a s . The latter results in a better fit 
but larger continuum value. 
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FIG. 10. Continuum extrapolation of 3 Pi — 3 Po- We use the same symbols as in Figure [)] to plot this element of the fine 
structure against the spatial lattice spacing determined from IP — IS. 
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FIG. 11. Higher excitation in bottomonium, T'. Here we plot the conveniently normalised splitting, 2S — IS, against the 
spatial lattice spacing. Within the errors we find a good agreement for all different methods and a slightly higher continuum 
limit than experiment. From potential models, we can expect that this ratio will be lower in dynamical simulations. 
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